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Optimum airfoils in the sense of maximum lift coefficient are obtained 
by a newly developed method. The maximum lift coefficient is achieved by 
requiring that the turbulent skin friction be zero in the pressure rise 
region on the upper surface. Under this constraint, the pressure dis- 
tribution is optimized. The optimum pressure distribution consists of 
a uniform stagnation pressure on the lower surface, a uniform minimum 
pressure on the upper surface immediately downstream of the front 
stagnation point followed by a Stratford zero skin friction pressure 
rise. When multiple-element airfoils are under consideration, this 
optimum pressure distribution appears on every element. The parameters 
used to specify the pressure distribution on each element are the 
Reynolds number Re^ and the normalized trailing edge velocity q^. 

The newly developed method of design computes the velocity distribution 

on a given airfoil and modifies the airfoil contour in a systematic manner 

until the desired velocity distribution is achieved. There are no 

limitations on how many elements the airfoil to be designed can have. 

Numerical examples of one- and two-element airfoils are given. The 

Cj^ values obtained range from 2 to 2.5. 
max 
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I. Iir.RODUCTlON 

Airfoils with high maximum lift coefficient (Cj^ ) are desirable for 

max 

certain types of aircraft. With their large , these airfoils can be 

max 

operated at a comparatively low speed while producing a large amount of 
lift force. This characteristic is of primary interest in the design of 
aircraft with short take-off and landing run capability. 

Since the value of is limited by the occurrence of boundary 

max 

layer separation on the upper surface of an airfoil where the fluid is 
subject to an adverse pressure gradient, the approach to obtaining a high- 
lift airfoil Ic to delay or entirely suppress the flow separation. In 
1920, Lachmann and Handley-Page designed the first high-lift wing which 
consisted of a conventional airfoil and a small additional element called 
a leading edge slat. This leading edge slat gives the fluid which passes 
through the gap between the slat and the main airfoil a high velocity. 
Consequently, the boundary layer which grows on the upper surface of the 
main airfoil has more momentum than it would have in the absence of the 
slat. This hlgh-momentum boundary layer can sustain a steeper adverse 
pressure gradient and hence delay the separation. Therefore, by placing 
the leading edge slat at an appropriate position, the airfoil can be 
operated at a higher angle of attack without flow separation and the lift 
force Is increased. An alternative way of obtaining more lift force is 
to put an additional element near the trailing edge of the airfoil. 

The so-called trailing edge flap does not delay the separation of boundary 
layer on the upper surface of the airfoil. It produces additional 
circulation and hence lift force by both Its presence and the velocity 
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It Induces on the main airfoil. Performance of airfoils which have leading 

edge slats and/or trailing edge flaps are available in the literature. 

With all these auxiliary elements to provide more lift force for a given 

airfoil, few attempts have been made to find the maximum lift coefficient 

which can be obtained on a single element airfoil. The first successful 

work in this area was accomplished by Liebeck and Ormsbee.^ They designed a 

series of airfoils porisessing large values of by requiring the 

max 

pressure on the upper surface of an airfoil to rise in such a way that the 
turbulent skin friction is zero whatever the pressure is increasing. These 
airfoils are optimum ir. the sense of maximum lift coefficient. Zero skin- 
friction means that the boundary layer is about to separate. Thus, any 
attempts to obtain a larger lift coefficient by either increasing peak 
velocity or shifting the point of peak velocity in the downstream direction 
will result in a flow s«;paration. With this zero skin-friction require- 
ment, a variational problem was set up in order to find the pressure 
distribution which would provide the maximum lift coefficient. The 
solution to this variational problem specified the pressure distribution 
on the upper surface of the airfoil and one of the standard methods of 
airfoil design was employed to obtain the geometry of the optimum airfoil. 

Since more lift can always be obtained by appropriately using a 
leading edge slat and/or a trailing flap, a natural extension of Liebeck 
and Ormsbee's work would be the search for optimum multiple-element 
airfoils in the sense of maximum lift coefficient. This is the goal of 
this research. The study consists of two parts. The first 
part deals with the search for pressure distributions which 
provide a maximum lift coefficient. The second part determines 
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the geometry of the multiple-element airfoil which produces theee preiiure 
distributions. As will be seen in the next chapter, the optimum pressure 
distribution is found to be identical to the one for single-element airfoils. 

A new method of two-dimensional airfoil design is developed to obtain the 
geometries of the airfoils which produce these optimum pressure distributions. 
Since no methods of imiltiple-element airfoil design are available in the 
literature, the method developed in this research may be considered at a 
major contribution to this aspect of aeronautical science. This method 
makes it possible to obtain the geometries of each element of the nniltiple- 
element airfoil when the pressure distribution is specified along the 
surface of each element. Numerical examples of two-element airfoils are 
treated using this new method of design. The fluid is considered to be 
Incompressible and Reynolds numbers are assumed to be sufficiently large 
that the boundary layer in the zero skin- friction state is turbulent. 
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II. DETERMINATION OF OPTIMUM PRESSURE DISTRIBUTION 
A. The Zero Skin-Friction Requirement and 

^max 

It is well Imown that the maxim>.3n lift coefficient attainable xor a 
given airfoil is determined by the flow separation %rhich occurs on the 
upper surface of the airfoil. Generally, the lift coefficient of an air* 
foil is directly proportional to its angle of attack with respect to the 
free stream. For angles of attack larger than some particular value, 
however, the lift coefficient increases at a slower rate with respect to 
the Increase of angle of attack and even starts to decrease when the angle 
of attack is increased beyond another particular value, <X2* reason 

for this slow increase of lift coefficient and a later decrease is that 
the boundary layer on the upper surface of the airfoil cannot sustain 
the steep pressure rise which appears on the upper surface as a result of 
a large angle of attack. At a flow separation appears at a point near 
the trailing edge and generates a small region of reverse flow. This 
region of reverse flow expands as the angle of attack becomes larger than 
and finally develops into a iarge high pressure region with a large 
portion of the airfoil upper surface as its boundary. 

While flow separation is a local phenomenon, it is apparent that when 
the flow starts to separate at a certain point on the airfoil surface, the 
boundary layers at other points may be still able to sustain their local 
pressure gradient. Therefore, from the viewpoint of generating large lift 
coefficient, it seems that an airfoil would be more efficient if it could 
generate a pressure distribution such that the boundary layer separation 
would occur ^.erever the pressure is increasing if separation does occur at 
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all. In other words, an airfoil which can provide a large maximum lift 
coefficient Is the one which has no flow separation for angles of attack 

smaller than a value a > and suddenly stalls at a ■ a • This point may be 

8 8 

made clear by considering the pressure distribution shown in Figure 1. 

The lift coefficient of the airfoil which genetates this pressure 
distribution Is represented by the area enclosed by the curve Cp(x) where 
the x-axls is parallel to the free stream. If the flow separates at A, 
but at no other point, then the skin friction has positive value for every 
point upstream of A. Now, let the pressure distribution between B and A 
follow a different i^urve such as the dotted line In Figure 1. Then, the 
pressure gradient at points upstream of A will be steeper and the area 
enclosed by the C^(x) curve is increased which means that the lift cO' 
efficient has been increased too. It is possible that the flow will 
separate somewhere between A and B V7lth this new pressure distribution. 

Let it be assumed that the skin friction resulted from pressure distribution 
® is zero everywhere between C and A. It is not hard to see that any 
attempts to increase the lift coefficient of the airfoil by forcing the 
pressure distribution to follow a curve such as the one represented by the 
broken line will steepen the pressure gradient and hence introduce flow 
separation. The lift coefficient, then, will be decreased. Therefore, 
among all the pressure distributions which are identical at every point 
except the region where pressure is increasing, the one which yields zero 
skin friction over the pressure recovery region provides the maximum lift 
coefficient. This is the foundation based upon which high lift airfoils 
of the Liebeck-Ormsbee type are designed. When an additional elment is 
present near the airfoil, such as a trailing edge flap, the argument 
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stated above still applies. Namely, the pressure recovery on both the 
main airfoil and the auxiliary element should yield zero skin friction in 
order to achieve maximum lift coefficient of the airfoil. In general, 
a multiple-element airfoil will have maximum lift coefficient if this 
pressure recovery of zero skin-friction type appears on each element of 
the airfoil. 

B. Variational Probl«n and Its Solution 

1. Stratford’s Pressure Distribution and Formulation of 
the Variational Problem 

As stated in Chapter I, the fluid flows are considered to be in- 
compressible. In order to obtain a boundary layer which is capable of 
sustaining a large positive pressure gradient, only flows with large 
Reynolds numbers are considered. Based on these two essential assumptions, 
a problem in the calculus of variations will be formulated in this 
subsection. 

Suppose that the fluid is inviscid. The lift force of an airfoil 
is given by the Kutta-Joukowski Theorem 

L - pn.r 

where p is the density of the fluid, is the speed at which the airfoil 
travels in the fluid and F is the circulation required to place the rear 
stagnation point at the trailing edge of the airfoil. With a given air- 
foil of specified chord length and a specified velocity, the only way to 
obtain more lift force is to increase the circulation F. The circulation 
F, by definition, can be written as 

®T 

F = J* v(s)ds 
0 


( 1 ) 
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where s Is the arc length of the airfoil contour, measured clockwise 
from the trailing edge along the airfoil surface (see Figure 2), s^ is the 
total length of the contour of the airfoil and v(s) is the surface 
velocity distribution. With subscript L denoting the front stagnation 
point of the airfoil, the Integral (1) can be decomposed into two terms 


L T 

J v(s)ds + J v(s)ds 
0 s. 


( 2 ) 


where 0 < s < s^ holds for the lower surface and < s < s^ holds for 
the upper surface of the airfoil. Notice that along the lower surface, 
the direction of v is opposite that of increasing s. Therefore, the first 
term in (2) is a negative one. There is no way that a positive quantity 
may be obtained as the first term of (2) unless the front stagnation point 
is forced to coincide with the trailing edge. However, this is physically 
unattainable. Thus, the most which can be obtained from the lower surface 
is a zero velocity all the way from the front stagnation point to the 
trailing edge. In other words, v(s) * 0 for 0 < s < Sj^. With this 
velocity on the lower surface of an airfoil, variations in lift may occur 
only by changing the velocity distribution on the upper surface. The 
circulation then stands as 

®T 

r ■ J v(s)ds . 

For convenience, the origin of s may be shifted from the trailing edge to 
the front stagnation point. With this modification, the circulation is 

®U 

J v(s)ds 
0 
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where and subscript U stands for the upper surface. The lift 

coefficient of the airfoil, based on free stream dynamic head and length 
Sy of the upper surface, is 

s 


or 


C 


L 


PIL 


2 '’"-“u 


u 

J v(s)ds 
0 


= f- J q(s)ds 
®U 0 


( 3 ) 


where q(s) is the surface velocity normalized with respect to the free 

stream velocity. This velocity q(s) always starts from the front 

stagnation point s = 0 with a value of zero, increases monotonically to 

a maximum value q^ at a point s » s^ and then decreases monotonically to 

the value of zero at the trailing edge where s = s^. For the purpose of 

generating maximum lift coefficient in the way described in Section A 

above, the function q(s) can be any form for 0 < s < Sq. There is no 

restriction for that part of q(s) other than the one that it must be 

monotonically non -decreasing. The function q(s) for < s < s^, however, 

must satisfy the requirement of being able to provide a boundary layer of 

zero skin friction everywhere between s^ and s^, with the reason described 

In Section A. Such a function q(s) may be obtained by modifying the 

2 

expression derived by Stratford. Stratford considered the turbulent 
boundary layer grown on a flat plate and derived an expression for pre- 
dicting the separation of this turbulent boundary layer when an adverse 
pressure gradient is encountered. Let the leading edge of the flat plate 
be at the origin x » 0 and let the pressure start to rise at x > x^. 

Separation will occur when the following equation is satisfied by C (x) 

0 




.X 
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(4) 


where C^ (x)ls the pressure coefficient based on free stream p^, U^; $ 
is a constant which is approximately 0.66 for Reynolds number of order 
10^; n is the C(xii3on logarithm of Re which is the Reynolds number based 
on Uq and x. The criterion employed in deriving (4) is that separation 
is imminent wherever local skin friction equals zero. Therefore, by treating 
(4) as a differential equation for C (x) , a pressure distribution which 


provides zero skin friction at every point downstream of x^ is found to be 

i i 

" = 0.645(0.435 [ (~)^ - l]}^^"for < -^ 


PO *0 


= 1 - 


(x+b) 


1/2 


for C > 

Po - n +1 


(5) 


where Re^ is the Reynolds number based on x^ and uniform velocity and 
n is the conanon logarithm of Re^; a and b are constants to be determined 


such that C and dC /dx are continuous at C 


n-2 

nfl* 


These derived 


results have been verified by Stratford and experiments showed that this 
C (x/x ) does provide a boundary layer of zero turbulent skin friction. 

Po ° 

The pressure starts to rise at s = s« on an airfoil where q = q_. 


Hence the relation between q(s) and (x/x^) is established by first 
setting Uq - Then, since 

U 


U U "0 . - 

— m — — and c 

U U„ U p„ 

09 0 * *^0 


- I 
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one obtains 


q(s) 


-oti 


'p <r'i 
^0 0 


1/2 


( 6 ) 


It should be pointed out, however, that expression (5) above was derived 
when the boundary layer on a flat plate in a uniform flow was under 
consideration. For a non-decreasing velocity distribution U(s) the 
relation between Sq and is 


X 


0 


J 

0 0 


ds 


(7) 


which merely states that the momentum thickness of the boundary layer at 
8 » Sq on an airfoil with velocity distribution U(s) for 0 < s < Sq has 
been set to be the same as that of the boundary layer at x * on a flat 
plate in a uniform flow. Now, with the small velocity near the front 
stagnation point, the boundary layer there is likely to be laminar. This 
laminar initial boundary layer is acceptable provided a transition to 
turbulent boundary layer has taken place before the pressure starts to 
Increase. This is because Stratford considered only the separation of 
turbulent boundary layers. With this laminar boundary layer present, 
equation (7) is superseded by 


38 . 2 (^) 

t t 


,Vl/8 

“t 


y 




’t 


V 


( 8 ) 


where subscript t indicates that the variables are evaluated at the 
transition point. An Instantaneous transition with the preservation of 
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moi&ent\]in thickness has been assumed in deriving equation (8). Therefore, 
a boundary layer which is entirely laminar between s « 0 and s * s^ is 
acceptable provided an instantaneous transition takes place at s^. In order 
to have pressure increase at s^, the starting point of q(s) will not be 
the same as the equivalent flat plate boundary layer. This is illustrated in 
Figure 3. Hence equation (6) must be used with the understanding that s 
and X are related by equation (8) and Figure 3. 

Now, with an arbitrary acceleration q(s) for 0 < s < s^ and a pressure 
rise of the Stratford's type described by (5), the lift coefficient of the 
airfoil is 


' 7" y - C (^)]^^^ds} 

^ ft L ^ Pq *0 


where is given by equation (5) with defined in equation (8) and the 


relation between x and s is shown i;i Figure 3. In this expression, s 


PO 


u 

can be considered a constant, s^ is not specified. q(s) is any mono* 
tonically increasing function, q^ Ls not specified and the function C 
depends on x/x^ and Re^ only. These quantities can be varied in order 
that a maxim\oa value of may be obtained. Among these quantities, Re^ 
should be specified independent of jhe others because Re^ specifies the 
boundary layer characteristics at s^^ and hence the zero skin friction 


pressure recovery (x/x^). When Re^ and s^ are given, the problem of 


finding the maximum value of becomes one of searching for a function 
q(s> and values of the quantities s^, q^ such that the quantity 




u 0 


0 '‘0 
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takes its maximum value, x^ Is to be computed from equation (8). This is 
a problem of calculus of variations for the funtion q(s). For q^ and s^, 
it is a problem of differential calculus. 

As can be easily seen from the results obtained by Stratford, the 
pressure recovery with zero skin friction has a large pressure gradient 
immediately downstream of x^. When the value of x is very large, the rate 
of Increase of pressure is very slow. In the limit, C will reach a value 
of one only If x goes to infinity. This means that the zero skin friction 
pressure recovery, when produced by an airfoil, will reach stagnation at the 
trailing edge only if the airfoil has an infinite chord length. In view of 
the fact that an airfoil with a chord of infinite length is not practical, 
a non-zero value of will be accepted as the trailing edge velocity. This 
leads to a sharp trailing edge which is acceptable for aerodyna.nic consider- 
ation. Now, the three quantities s^, q^ and q^ are related to each other 

once an Re and an s,. are given. Only the value of one of thm may be 
o u 

specified and the solution to the variational problem will determine the 
others. Suppose the value of s^ is specified. Then q^ can always be 
obtained for a given q^. The consequence is that there are infinitely many 
sets of q^ and q^ which can yield zero skin friction pressure recovery. 

The maximum lift coefficient which may be obtained under this circumstance 
has no upper bound. This indicates that the problem is not well defined. 
Alternatively, if q^ is specified, the lift coefficient will take a 
maximum value only when s^ goes to zero, and this solution does not make 
any physical sense. Therefore, the only remaining possibility is to 
specify the value of q^ and let the values of q^ and s^ be determined by 
the solution of the variational problem. Therefore, the problem can be 
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redefined as one of finding a function q(s) for 0 < s < Sq and two 
quantities q^ and s^ In order that 



has a maximum value. The values of Re^, s^j and q^ are specified and 
Is to be computed from equation (8). Since equation (8) takes a simpler 
form If the boundary layer Is turbulent from s > 0 to s ■ s^, the case of 
having a laminar Initial boundary layer will be considered after the 
simpler case of an all turbulent boundary layer Is treated. 


2. Solution to the Variational Problem When Boundary Lave r 
Is All Turbulent 


When the boundary layer Is all turbulent, equation (8) becomes 


X 


0 



ds 


( 10 ) 


For the convenience of algebraic manipulation, it Is desirable to make a 
transformation on the variable s. With the definition of a constant 
K > relation between s and x is 

s-x+(k-l)xQ 

or 

X - s - (k - 1) Xq . 

Dividing X by Xq, a variable z may be defined as 

z - ^ 

^0 

- J- - (k - 1) . 

0 
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with this notation, equation <9) may be written as 


" z+fcr ^ J ‘*^*^*** ‘‘O J <**3 (11) 

-(k-1) 1 0 


where Z denotes the value of z at the trailing edge. Now, the function C 


may take only one form of equation (5) or both depending on the values of 
j- and x/Xrt. For Re„ 

u u u p. 


Re^ and x/x^. For Re„ of order 10^ - 10®, will reach the value 


for a fairly small value of x/x.. Hence it can be assumed that C will 

n-2 ’■O 

take the value at a point s^ which is between Sq and s^ for the 

Reynolds numbers considered in this research. With this assumption, 

2 

Stratford's equation (29), which expresses the function C (x/x^) for 

n-2 ’’® 

C > ~rr » can be modified to give 
Pq - n+1 


‘0 '*U 


Sim 


1/4 


/a’ 


( 12 ) 


where a' ■ a/Ag? b* > Hence equation (11) becomes 

2 i .... 2 

“ Z+k-1 J 9 (*)<** + 

1-k /a 


m 


X J* {1-0.645C0.435 ReQ^^^(z^^^-l)]^^"}^^^ dz 
1 


q„(Z+b')^^'^ J* 

z+*t 1 u ^ (z+b') 

n 


dz 


(13) 


where Z ■> s /x. - (k-1), and equation (10) beccxnes 
m m u 


J dz . 

1-k ® 


( 14 ) 
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To find a stationary value o£ C^, the variation of C^ caused by the 
variation of the right hand side of equation (13) must vanish. This 
variation includes the variation of q(z), the variation of k and the 
variation of Z. The variation of k is related to that of q(z) by the 
constraint equation (14). Hence only the variation of q(z) and Z need to 
be considered. The appearance of q(z) inside an integral sign makes the 
problem fall into the category of calculus of variations. However, due to 
the simplicity of the integrands which contain q(z) in both equation (13) 
and equation (14) , the solution to this problem of calculus of variations 
can be obtained in a quite simple way. Because no derivatives of q(z) 
appear in the integrand, the Euler's equation resulted from the vanishing of 
first variation of is an algebraic equation instead of a differential 
equation. The solution to that equation is merely q(z) = constant. This 
constant, then, is determined by the requirement' that q(sQ) = qQ* Sub- 
stitution of q(z) = q^ into equation (14) gives k = 1. Therefore, as far 
as q(z) is concerned, Cj^ will have a stationary value only when q(z) = q^. 
The fact that the value of provided by this q(z) is really a maximum 
one may be established in the following way. Suppose a permissible 
variation is introduced into this function q(z). Permissible means that 

q(z) is still a monotonically non-decreasing function. The value of k 

1 

will always increase. The value of J q(z)dz may increase or it may 

1-k 

decrease . Because of the constraint of equation (14) , the increase of k 
is so large that C. will always be decreased by this variation of q(z). 
Therefore, in order to obtain a maximum lift coefficient, the fluid 
should be accelerated abruptly from the front stagnation point to a 
velocity qQ and remain at that value until pressure starts to increase. 
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With this conclusion, equation (13) becomes 


2 qy(Z+M^^* 2 


L Z 


/a' 


/a’ 


X J il-0.645[0.435 dz 


2 2 
■^2 ^ . J 


/a* Z_ (z+b’)^^^ 

ffl 


dz . 


Since depends on Re^ only, the integral in the second term may be 
abbreviated as 

KRCq) = J {1-0.645[0,435 dz . 

1 

This allows to be expressed as 

2q„(Z+b')^^^ . 

C = _U (1 + i(Re ) + I /a'[(Z+b’)^^^-(Z +b')^^^]} 

where the integration appearing in the last term has been carried out. 
With values of q^ and Re^ given, this varies with Z only , The maximum 
value of may be obtained by evaluating at a value of Z where the 
first derivative of C, with respect to Z vanishes. Taking the first 
derivative of C with respect to Z and setting the result equal to zero 
gives 

I [1 + KRejj) - I /a' (Z^+b')^''^l(Z+b'> + | /a'b’ (Z+b')^^^ 

+ ^ [1 + KRe^) - I /a* (Z^+b’)^^^] - 0 . (15) 
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1/4 

This is a fourth degree algebraic equation for (Z+b') and there are 

2 2 

four roots. The second derivative d must be evaluated at each root 

of equation (15) in order to find out which root gives the maximum value 

2 2 

of Cj^. The expression for d C^^/dZ is 


/a' 


{[1 + 


I(Reo) 


-f/.' 


(Z +b')- 
m 


3/4^^ 2(Z+b') 


1/4 


(Z±^ 


2Z 


- 3/4 


3(Z-<-b') 

16Z 



+ 


8 

3 


/a’b' 



(16) 


3 8 

For 5 X 10“’ < Re^ < 1 X 10 , equation (15) has two conjugate complex roots, 

one negative real root and one positive real root. In view of equation (12) 

1/4 

only positive values of (Z+b') can be accepted. CompUwations have shown 

that this positive root of equation (15) gives a value Z which is larger 

2 2 

than Z and this Z does yield a negative value of d C. /dZ , Substituting 
m t 

this Z into equation (12) gives the value of Qq/^I^ ®"d thus the solution 

to the variational problem is complete for an all turbulent boundary layer. 

(See Figure 4). The dependence of Z^, Z, and qQ/^u is given in 

2 2 

Table 1 . Also shown in Table 1 are the values of d C^/dZ . 

3. Solution to the Variational Problem When the Initial Boundary 
Laver Is Laminar 

When a laminar initial boundary layer is present, the full expression 
(8) must be used as the constraint on q(s) for 0 < s < Sq. Using the 
definition of q, equation (8) can be rewritten as 



( 17 ) 













Table 1 


Dependence of Z^, Z, ^ ^®0 ***' 

Turbulent Boundary Layer 


R6q '/ 10"^ 

Z 

m 

Z 

^0^% 

2 2 
d% /dZ^ 

5 

1.6625 

4.2340 

2.0822 

-0.003332 

6 

1.6502 

4.0713 

2.0824 

-0.003887 

7 

1.6400 

3.9473 

2.0835 

-0.004395 

8 

1.6313 

3.8487 

2.0852 

-0.004863 

9 

1.6236 

3.7677 

2.0872 

-0.005298 

10 

1.6168 

3.6995 

2.0893 

-0.005705 

20 

1.5735 

3.3312 

2.1099 

-0.008790 

30 

1.5493 

3.1640 

2.1264 

-0.010926 

40 

1.5327 

3.0611 

2.1397 

-0.012587 

50 

1.5201 

2.9887 

2.1508 

-0,013959 

60 

1.5101 

2.9338 

2.1603 

-0.015134 

70 

1.5017 

2.8900 

2.1686 

-0.016164 

80 

1.4945 

2.8538 

2.1761 

-0.017085 

90 

1.4883 

2.8233 

2.1828 

-0.017917 

100 

1.4827 

2.7969 

2.1889 

-0.018678 

200 

1.4478 

2.6430 

2.2318 

-0.024071 

300 

1.4285 

2.5662 

2.2587 

-0.027529 

400 

1.4154 

2.5165 

2.2784 

*0.030114 

500 

1.4054 

2.4804 

2 . 2942 

-0.032193 

600 

1.3974 

2. ^^523 

2.3073 

-0.033938 

700 

1 . 3908 

2.4295 

2.3185 

-0.035447 

800 

1.3852 

2.4104 

2.3283 

-0,036779 

900 

1.3802 

2.3940 

2.3371 

-0.0J7970 

1000 

1.3759 

2.3797 

2.3450 

-0.039051 
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Now, by defining g > replacing s by z, equation (17) becomes 

1 . 38.2 t''r" 

x«u q^ '■ J q« •' 




1-k 


■ f dz 

1-k-fgk 


(18) 


Since the expression for C^, equation (9), holds true whether the Initial 
boundary layer is laminar or turbulent, the problem is again to find a 
q(z) which will give a maximum value of 

1 


C, = 


L Z+k 


7 J ‘•(^>■*"■"55^: 


2 q„(2+b’)‘^* 


1-k 


/a' 


-m 


X J {l-0.645[0.435 ReQ^^^(z^^^-l) dz 
1 




dz 


Z (z+b')^^^ 
m 


(13) 


It can be seen that no derivatives of the function q(z) appear in the 
integral terms in equation (13) and equation (18). Thus, the solution 
to this variational problem is again q(z) = q^ for 0 < s < s^. With this 
q(z), equation (18) gives 


1 - 38. 2(—^)^^® (gk)^^®+ (l-g)k . 


(19) 


The critical Reynolds number at which boundary layer transition takes 
place is defined by 


Re 


»t"o 


cr 


V 
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r .i- -^--.*f1P.^^,S^..#r-'‘ 


which is 


Re 


cr 


8“o“o 


«'“o“o 


Hence , 


gk 


Re 

t 

Re, 


cr 


Substituting this gk into equation (19) gives 


1-38.2 Re."^ Re + -=r^ 
0 cr Re,, 


( 20 ) 


This result requires that a critical Reynolds number Re^^ be specified 
in addition to the Reynolds number Re^ when the initial boundary layer 
is laminar. Notice that the case Re > Re,, need not be excluded because 
an all-laminar boundary layer between s = 0 and Sq will give k > 1 and 
g « 1. When the boundary layer is all turbulent g ■ 0 and equation (19) 
gives k - 1 which is the result obtained in the previous subsection. 

With q(z) = q^ and k given by (20), equation (13) can be written as 

- drr Out'' + 1 


/a' 


X [(Z+b')^^^ - (Z +b')^^^]} (21) 

nri 


where Z, I(Re^) and all have the same meaning as In the previous sub- 
section. When the laminar initial boundary layer is absent, g > 0 gives 
k » 1 which reduces equation (21) to the form assumed by an all turbulent 
boundary layer given in Subsection 2. This verifies the consistency of 
the expressions derived up to this stage . 
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To find the maximum value of C with various Z values, the derivative 

L 

2 2 

dC, /dz is set to zero and the second derivative d C /dZ must be evaluated 
L Li 

at certain Z values. The equation which results from setting dCj^/dZ to 
zero is 

7 [k + l(Re^) - I /a' (Z +b' )^^^](Z+b' )- ^ /a' (k-l-b’ ) (Z+b’ ) 

4 U J m j 

- i (k-l-b')[k+l(Rey) - I /a' (Z^+b')^^^l * 0 . (22) 

As in the equation obtained in Subsection 2, this is also a fourth degree 

algebraic equation for the variable (Z+b')^^^. For 2 X 10^ < Re^ < 1 x 10^, 

two different values have been used for Re . They are 5 x 10^ and 

cr 

1 X 10^. In both cases, there are always two conjugate complex roots, 

one negative real root and oi e positive real root. As stated in the 

previous subsection, only i;h>-. positive real root is retained, and it does 

2 2 

jive a negative value of d when this positive real root is sub- 

stituted into 


4 

— {[k + KRe^,) - /a' 

/a' ° ^ 


(Z +b')^^^]C- 
m ^ 


3_ (Z+b*)~ 
16 k+Z-1 


JM^Ul2!1+ 2(Z-:-b»)^^^ ^ 

2(k+Z-l)^ (k+Z-1)^ 


8 /a' ^ G /a* (Z+b* ) ■> 

*1 2 ' 3 ^ 

^ (k+Z-i) ' (k+Z-l)^ 


This means that a maximum value of does exist at this Z. Substituting 

this Z into equation (12) gives the value of qg/q^* The dependence of 

0 . Z . Z and q»/q.. on Re^ are shown in Table 2 and Table 3. Also shown 
* ’ uj ^0 , U 0 
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Table 2 

2 2 

Dependence of g, k, 2 , Z, qo/q„ and d C, /dZ on Re_ for 

m u u i-i V ^ 

Boundary Layers Which are Initially Laminar, = 5 X 10 


0 >' I”''* 

g 

k 

Z 

m 

Z 

’obu 

2 2 
d^C^/dZ 

15 

0.979 

3.404 

1 . 7462 

6.7252 

2.1895 

-0.002548 

20 

0.892 

2.804 

1.7258 

6.0092 

2.1656 

-0.003025 

25 

0.819 

2.443 

1.7101 

5.5548 

2.1512 

-0,00341 7 

30 

0.757 

2.202 

1.6974 

5,2363 

2.1419 

-0.003754 

35 

0.704 

2.031 

1.6868 

4.9983 

2.1355 

-0.004053 

40 

0.657 

1.902 

1.6776 

4.8121 

2.1331 

-0.004323 

45 

0.617 

1.802 

1.6696 

4.6613 

2.1280 

-0.004573 

50 

0.581 

1.721 

1.6625 

4.5361 

2.1257 

-0.004806 

55 

0.549 

1.656 

1.6561 

4,4299 

2.1241 

-0.005026 

60 

0.520 

1.601 

1.6502 

4.3384 

2.1230 

-0.005235 

65 

0.495 

1.555 

1.6449 

4 . 2584 

2.1222 

-0.005435 

70 

0.471 

1 .515 

1.6400 

4.1876 

2.1217 

-0.005627 

75 

0.450 

1.481 

1.6355 

4.1245 

2.1214 

-0.005811 

80 

0.431 

1 .451 

1.6313 

4.0677 

2.1213 

-0.005990 

85 

0.413 

1.424 

1.6273 

4.0161 

2.1213 

-0.006163 

90 

0.397 

1.401 

1.623b 

3 96^1 

2.1215 

-0.006331 

95 

0.381 

1.380 

1.6201 

3.9260 

2.1217 

-0.006494 

100 

0.367 

1.361 

1.6168 

3.8863 

2.1219 

-0,006653 

200 

0.212 

1.180 

1.5735 

3.4410 

2.1236 

-0.009236 

300 

0.149 

1.120 

1,5493 

3.2423 

2.1441 

0.011152 

400 

0.115 

1.090 

1.532? 

3.1219 

2,1543 

-0.012697 

500 

0.093 

1.072 

1,5201 

3.0385 

2,1633 

-0,013998 

600 

0,079 

1.060 

1,5101 

2,9759 

2,1712 

-0.015127 

700 

0.068 

1.052 

1,5017 

2,9264 

2.1784 

■0.016127 

800 

0.060 

1.045 

1.4945 

2.8860 

2.1849 

-0.017026 

900 

0.053 

1.040 

1.4883 

2.8520 

2.1909 

-0.01 7842 

1000 

0.048 

1.036 

1,4827 

2.8229 

2.1964 

-0.018592 
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Table 3 

2 2 

Dependence of g, k, Z^, Z, <1 q/<1u ^®0 

Boundary Layers Which are Initially Laminar, * 10^ 


RBq X lo"^ 

g 

k 

Z 

m 

Z 

lo'iu 

2 2 
d Cj^/dZ^ 

25 

0.966 

4,141 

1.7101 

6.8501 

2.2816 

•0.002902 

30 

0.922 

3.617 

1.6974 

6.3104 

2.2585 

-0.003317 

35 

0.881 

3.243 

1.6868 

5.9131 

2.2412 

-0.003689 

40 

0.843 

2.963 

1.6776 

5.6070 

2.2277 

-0.004026 

45 

0.810 

2.745 

1.6696 

5 . 3630 

2.2171 

-0.004336 

50 

0.778 

2.570 

1.6625 

5.1634 

2.2085 

-0.004621 

55 

0.749 

2.428 

1.6561 

4.9965 

2.2014 

-0.004886 

60 

0,722 

2.309 

1.6502 

4.8547 

2.1955 

-0.005133 

65 

0.697 

2.208 

1 . 6449 

4.7323 

2.1906 

-0.005366 

70 

0.673 

2.122 

1.6400 

4.6255 

2.1865 

-0.005587 

75 

0.651 

2.047 

1.6355 

4.5313 

2.1829 

-0.005796 

80 

0.631 

1.981 

1.6313 

4.4475 

2.1799 

-0.005995 

85 

0.612 

1.924 

1.6273 

4.3723 

2.1773 

-0.006185 

90 

0.593 

1.872 

1.6236 

4 . 3044 

2.1751 

-0.006368 

95 

0.577 

1.827 

1.6201 

4.2427 

2.1732 

-0,006544 

100 

0.560 

1.785 

1.6168 

4.1864 

2.1715 

-0.006713 

200 

0.359 

1.393 

1.5735 

3.5887 

2,1621 

-0.009317 

300 

0.264 

1.262 

1.5493 

3.3415 

2.1660 

-0.011179 

400 

0.209 

1.196 

1.5327 

3.1970 

2.1719 

-0.012677 

500 

0.173 

1.157 

1.5201 

3.0990 

2.1782 

-0.013944 

600 

0.147 

1.131 

1.5101 

3,0266 

2.1842 

-0.015049 

700 

0.128 

1.112 

1.5017 

7.. 9702 

2,1900 

-0.016030 

800 

0.114 

1.098 

1.4945 

^ 9244 

2.1953 

-0.016916 

900 

0.102 

1.087 

1.4883 

2.8863 

2.2004 

-0.017722 

1000 

0.092 

1.079 

1.4827 

2.8538 

2.2051 

-0.018465 
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2 2 

are the valuc-3 of d C^^/dZ . The velocity distribution q(s' for 
0 < 8 < is sketched in Figure 5. Notice that the location where 
pressure starts to increase is different from the one shown in Figure '■* 
This is due to the presence of the laminar Initial boundary layer 


C. The Pressure Distribution To Be Used in Designing Two-El<-m e nt 
Optimum Airfoils 

The optimum pressure distribution which will provide ar airfoil with 

maximum lift coefficient, as shown in Section 6, consists of a const a*'*' 

stagnation pressure along the lower surface of the airfoil, and an atrup' 

pressure drop to a minimum value which covers a certain distance on thr 

upper surface of the airfoil until the pressure starts to increase 

according to Stratford's zero skin friction pressure distribution, The 

determination of geometries of multiple-element airfoils is the main 

object of this research. Pressure distributions will now be considered 

which are optimum with respect to for this case. The lift force 

max 

is again pU T according to the Kut ta-Joukowski Theorem where T is the 
total circulation of the airfoil and can be written as 



r = J v^fSjjdSj + J v^(s^>ds^ + J V2(s2.)ds2 



V 2 ( 82 .'>ds 


2 


for airfoils consisting of two elements. Contour length is again 
denoted by s which Is measured clockwise from the trailing edge of -ach 
element; subscript 1 denotes element one and subscript 2 denotes element 


imized velocity distribution on the upper 
face of a single-element airfoil with 
linar initial boundary layer. 





tv70. The fiiir: -.anv. aaa th<; third term are n^igative quantit :t33 imce 
the surface velocity is opposite the direction of increasing s. The 


maximum value of these quantities is zero. Thus, only two terms are left 
ana the lift coefficient may be written as 

h'ii' ■^1 J' '»2'h>'‘h 


where i is a characteristic length. For the convenience of manipulation, 
s„ 




s can be taken as Jt. This is not only convenient, it is also a 

L* > 


moce-or-lass crcioditional w-ay because s - s is appr*jximately the chord 

1 h 

length of element one. Wnon the origin of measuring s is moved to the 
front stagnaticr. point of each element as was dons in Subsection B.l, 


C becomes 

la 


= 


- 1 

^1 qiC-iMs, + -- 

0 ''i ‘•'z 0 


J ■ 


The first term represents tne contribution of element one to , and it 
has exactly the same form as equation (3). The second term represents 
the contribution of element two to , and it has the same form as 

li 

equation (3) except for a multiplicative factor s /s,. . C is to be 

2 ^ 

maximized under the condition that pressure rise on upper surface of 
each element must yield a boundary layer of zero skin friction in the 
region wnero pressure is _ncr3asing. Th:.t: requiiement of obtaining 
iiiaximum l:.f': coefficient has been established in Section A. The 
pressure rise of zero skin friction is specified by a Reynolds number 
ke^. In the case of two-element airfoils, the specification of Re^ on 
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one element is entirely independent of that on the other. Furthermore, 
the pressure distribution on element one can be specified without con- 
sidering the pressure distribution on element two. la other words, the 
pressure distributions of two elements are independent. This fact can be 
easily seen when one realizes that the pressure distribution is specified 
along the surface contour and the relative ^/osition of the two elemerts 

is a free parameter at this stage. Therefore, assuming s, /s to be a 

u * u , 


constant, C will have a maximum value when both 


u. 




and 


'-J 

1 0 

S 




U, 


^ ^2 0 


take their maximum values. With Re^ and given on each element, the 
problem of maximizing Cj^ is equivalent to the problem of maximizing the 
lift coefficients of two single-element airfoils. The solution to this 
problem has been presented in Subsections B.2 and B.3. This means that 
for an air foil consisting of two elements, the maximum lift coefficient 
will be obtained when the pressure distribution on the surface of each 
element is such that it is an optimum pressure distribution if each 
element were considered to be a single-element airfoil with same values 
of RBq and q^ (see Figure 6). 

This result can be extended to airfoils consisting of any number of 
elements. Namely, for a multiple-element airfoil, the maximum lift 
coefficient will be obtained when the pressure distribution on each 
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element is composed of a constant stagnation pressure along the lover 
surface, an abrupt drop to a minimum pressure at the front stagnation 
point followed by a constant minimum pressure which covets a certain 
distance on the upper surface until the pressure starts to rise according 
to Stratford's zero skin friction pressure distribution. It can be re- 
called that the characteristic length s^ for a single element airfoil 
will not be available until values of p, U^, Re^ are specified. Jence 
s^ is part of the solution to the variational problem. This confirms 

the validity of assuming s, /s,, to be a constant when maximizing of 

U2 L 

a two-element airfoil. Table 4 shows the dependence of s^^ on Re^ for 

“6 2 

l^= 200 fps., p,/p = 160 X 10 ft /sec. This makes it possible to 
design two-element airfoils i^ose elements have the desirable chord 
ratios. 

Although the pressure distribution sho«m in Figure 6 provides an 
airfoil with a maximum lift coefficient, some modifications have to be 
made in order to realize a physically meaningful airfoil. First, the 
constant stagnation pressure along the lower surface is not possible to 
obtain. Since the lower surface of an airfoil does not contribute much 
to the lift force, any modification of the constant stagnation pressure 
distribution will not reduce the C value too much from its original 

la 

value. However, a pressure distribution which is monotonically de- 
creasing from leading edge to trailing edge is preferred, and a linear 
form for q is chosen in this study. A sionotonically decreasing pressure 
from leading edge to trailing edge will produce a boundary layer which 
always attaches to the airfoil surface, and a linear relation q(s) Is 
chosen for its simplicity. As will be seen In the next chapter, the 
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Table 4 

, on Re. for U ■ 200 £ps. , 

0 *2 
p ■ 160 X 10“6 ft^/sec. 

Boundary Layers Which are Initially Laminar 


Re^xlO*^ 

Re^^ - 5x10" 

CT 

8u(ft.) 

R.„ - 10 

15 

0.5004 


20 

0.5772 


25 

0.6506 

0.8758 

30 

0,7215 

0.9487 

35 

0.7905 

1,0190 

40 

0.8580 

1.0874 

45 

0.9242 

1.1542 

. 50 

0.9893 

1.2196 

55 

1.0535 

1 2840 

60 

1.1168 

1.3475 

65 

1.1794 

1,4101 

70 

1.2413 

1.4720 

75 

1.3025 

1.5332 

80 

1.3632 

1 5939 

85 

1.4234 

. 1.6540 

90 

1.4831 

1.7137 

95 

1.5424 

1.7728 

100 

1.6012 

1.8316 

200 

2.7169 

2.9463 

300 

3.7638 

3.9926 

400 

4.7713 

4.9995 

500 

5.7517 

5.9794 

600 

6.7117 

6.9388 

700 

7,6554 

7,8820 

800 

8.5856 

8.8116 

900 

9.5045 

9.7299 

1000 

10.4134 

10.6383 
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pr«s8ur« distribution along the airfoil contour Is specified In terms of 

q, not C . Therefore, the modified pressure distribution Is a linear 
P 

q(s) from the front stagnation point where q ■ 0, to the trailing edge 
where q ■ q^> However, the requirement that the airfoil produce such a 
linear q(s) Is not very crucial because the pressure distribution on the 
upper surface dominates the problem, and specification of q(s) over the 
entire airfoil contour may result in a geometry which Is not physically 
meaningful. For example, the upper surface might cross over the lower 
surface at a point between the leading edge and the trailing edge. Thus, 
during the course of determining the geometry of the airfoil, liberty is 
taken with the pressure distribution on the lower surface in order to obtain 
a physically meaningful gemaetry. In regard to the pressure distribution 
on the upper surface, a slight modification also must be made. This 
modification Is to change the abrupt pressure drop at the front stagnation 
point to a gradual one. Because an abrupt pressure drop at the front 
stagnation point corresponds to a leading edge of zero radius of curvature 
which Is not suitable for operating the airfoil at various angles of 
attack, this modification Is necessary. Knowledge about various types 
of airfoils will help in making a decision as to what is the most 
adop table form of this modification. 
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III. DETERMINATION OF THE GEOMETRIES OF AIRFOILS WHICH 
PRODUCE THE OPTIMUM PRESSURE DISTRIBUTION 

A. Review of Methods of Airfoil Design 
1. Sato's Method and Weber *s Method 

The methods of designing airfoils with velocity specified along the 
* 

airfoil surface fall into two categories. One is the method of conformal 
transformation and the other is the method of distribution of 
singularities. The most powerful method which belongs to the first 
category is the one developed by Sato. The basic formulae used by Sato 
actually are ::hose developed by Lighthill,^ but Lighthill's method has 
the disadvantage that velocity must be specified in terms of closed form 
functions in order to be able to carry out the necessary integration. 

With the help of high speed computers, Sato's method allows a velocity 
distribution of any kind to be specified, and the integrations are carried 
out numerically. The expression for the velocity distribution is assumed 
in such a way that the front and the rear stagnation points can be treated 
separately, A well-behaved function g(6) takes up the velocity distribution 
everywhere with the exception of the stagnation points and three constants 
A, B and T which are embedded. These constants are determined by the 
function g(@), the closure condition of the airfoil and the fact that flow 
field at infinity is a uniform one. Consequently, the resulting airfoil 
is always a closed curve and the disturbance dies out at large distances 
from the airfoil. The results are very accurate near the leading and the 
trailing edges where most other methods have difficulties. This is 
because the singularities caused by stagnation points are treated 
separately in an anaXv'tlcal way. The method is to be used on an 
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iterative basis in the sense that the computation must be repeated many 
times before a satisfactory geometry of the airfoil may be obtained. 

The reason of doing this is as follows A set of initial values must be 
given to A, B and t in order to obtain the function g(9) from the specified 
velocity distribution. This g(6) is to be used to obtain another set of 
A, B and t which will give a t.’oso'^. curve as the airfoil geometry. Unless 
the initial set of A, B, and t happens to hit the solution of the design 
problem, these two sets of A, B and t will not be the same and the re- 
sulting airfoil will not produce the desired velocity distribution. The 
function g(0), then, is modified in such a way that the airfoil geometry 
obtained in the next cycle will produce a velocity distribution which is 
closer to the desired one. This procedure is repeated until the specified 
velocity distribution is reached Therefore, Sato's method always 
guarantees an airfoil which produces the desired velocity distribution 
to be obtained. 

Turning to the method of distribution of singularities, one observes 
that the crudest one is the inverse of thin airfoil theory. Since large 
disturbances are not allowed in thin airfoil theory, airfoils with large 
thickness/ chord ratios and/or large cambers which are operating at high 
angles of attack can not be obtained by this method. Modifications have 
been made by Weber^ to include this capability by considering both first 
order and second order terms . 

2. Method of Conformal Transformation and Method of Distribution 
of Singularities in Designing Ifctlt ip le- Element Airfoils 

Vfhen an airfoil of more than one element is to be designed, both 
Sato's and Weber's methods fail. Because Sato's method en^loys conformal 
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tr«ns formation, the reason why it fails for multiple-element airfoils is 
obvious. Take an airfoil of two elements as an example. The domain 
outside this airfoil is a triply connected one while the domain outside 
a single>element airfoil is a doubly connected one. Therefore, all the 
theories and formulae developed for the latter cease to be applicable 
to the former. Weber's method seems more likely to be applicable to 
multiple-element airfoil design. But, using a conventional airfoil with 
a slotted flap as an example, the angle of attack and camber of the flap 
are not properly defined and their contributions to the velocity dis- 
tribution are difficult to identify. In view of the fact that no other 
methods of airfoil design are capable of treating airfoils consisting of 
more than one element, the necessity of developing a new method becomes 
clear. 

The first method considered was the conformal transformation. The 
reason for proceeding in this direction was that Sato's method shows that 
single-element airfoils with satisfactory geometry which produce the 
desired velocity distributions can always be obtained by employing a 
conformal transformation and modifying the transformation function in a 
systematical way. As shown by Garrick,^ the domain outside two closed 
contours can be transformed into the annular region between two concentric 
circles. Alternatively, a domain of rectangular shape may be obtained 
with the help of a logarithmic function. In this case, the two contours 
are mapped into two sides of the rectangle i^ich are facing each other. 
Garrick performed this transformation on two HACA 4412 airfoils and 
computed the velocity diatrlbutlon on the surface of each airfoil in 
the same fashion as Theodorsen computed the velocity distribution on the 
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surface of a single airfoil by transforming tue airfoil into a cricle. 
Therefore, by inverting Garrick's method, geometries of two-element air- 
foils should be obtainable by specifying the velocity distribution along 
the boundaries of the annular region between the concentric circles. 

Based upon this concept, a set of formulae, which states the closure 
condition of eabh element of the two>element airfoil, was derived. 
Expressions for velocity distributions were the same as the one used by 
Sato. It was impossible to obtain the explicit formulae for computing 
the constants. Instead, all the six constants which are part of the 
velocity distribution and hence part of the transformation function appear 
implicitly in six integral relations. When the numerical calculations 
were carried out on a computer, it was found that a tremendous amount of 
time and work was required to find one set of constants. Since this 
method is to be used also in an iterative way, the time and work involved 
in computing several sets of constants make this treatment of the design 
problem formidable. As indicated by Garrick, the relative position of 
the two elements varies with the values of four other constants, and no 
study has been made as to how these values determine the relative position. 
Because of these two drawbacks, the method of conformal transformation 
was considered intractable in designing multiple-element airfoils. Next 
to be considered was the method of distribution of singularities. 

The method of distribution of singularities was originally a method 
of computing the pressure distribution on the surface of a given airfoil. 
Singularities of unknown strengths such as sources, sinks, doublets or 
vortices are distributed inside the airfoil contour or on the airfoil 
surface. The strengths of these singularities are obtained by conq^uting 
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the velocities induced by these singularities end requiring the tangcncy 
condition be satisfied. The surface pressure distribution is then 
confuted from the strengths of these singularities. When this method 
is applied to solving a design problem, an iterative procedure is 
necessary. The starting point of this iterative procedure is an airfoil 
with arbitrary shape. The pressure distribution on the surface of this 
starting configuration is computed by the method of distribution of 
singularities. Then a modification of the geometry is performed according 
to how much the computed pressure distribution differs from the desired 
one. Generally, more than one modification is necessary and the iterative 
process continues until a satisfactory pressure distribution is achieved 
by the airfoil. At this stage, there are two questions which must be 
answered. One is how the pressure distribution on the surface of a given 
airfoil should be computed. The other is how to modify the geometry of 
the airfoil in such a way that the velocity distribution on the surface 
of the modified airfoil will be closer to the desired one than the 
previous one is. The answers to these two questions are presented in 
the following sections. 

B. Methods of Computing Pressure Distribution on the Surface 
of a Given Airfoil 

There are many methods of distribution of singularities available in 
the literature which coaq>ute the velocity distribution on the surface of 
a given airfoil. In the earlier methods, singularities are distributed 
inside the airfoil contour. The disturbances produced by the airfoil 
are considered to be composed of those due to thickness of the airfoil 
and those due to c^>er and «igle of attack. Sources and sinks are used 
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to represent the former, and the latter are represented by vortices. In 
the newly developed methods, all the singularities are distributed on 
the airfoil surface, and the pressure distributions are computed by 
solving integral equations of different forms. The outstanding feature 
of these methods is that the airfoils under consideration can possess any 
geometry and any orientation relative to the free stream. Effects of 
thickness, cambpr anl angle of attack do not have to be considered 
separately. As s. general property and hence a limitation of the method 
of distribution of singularities, the fluid flow is considered to be 
incompressible and nonviscous. 

Among all the new-fashioned methods, perhaps the best known is 

g 

the method of Hess and Smith. Although this method has gained such 
publicity that it almost becomes a standard method of computing velocity 
distribution on the surface of an airfoil in incompressible potential 
fluid flow, it has some drawbacks which have not been noticed by many 
people. There are two other methods which are very affective but not 
known to many people in this country. They were developed in Germany 
and one is by Martensen^ and Jacob, the other is by Oellers.^^ These 
three method.' will be analyzed and compared with each ot''^r in the 
following subsections. 

1. Hess-Smlth*s Method 

In this method, the airfoil surface is replaced by a source sheet 
with strength 9(s), where s is the distance measured along the airfoil 
surface. Considering the airfoil to be stationary, velocities induced by 
the source sheet are combined with the free strema velocity to satisfy the 



tangency condition that there should be no velocity component normal to 
the airfoil surface. This concept results in a Fredholm integral 
equation of the second kind 

2^(s) + J ff (s') ^ r(s,s')ds' » F(s) (23) 

0 

where denotes the trailing edge of the airfoil, r(s,s') is the distance 
between two poijits represented by s and s', and F(s) is a function related 
to onset flows. For a uniform free stream, F(s) * - * n(s) where n(s) 

is the local outward unit normal vector. In this equation, the first 
term on the left hand side represents the normal velocity at point s 
induced by the local source 0(s). The second term represents the normal 
velocity at point s Induced by the remaining source sheet. When the 
geometry of an airfoil is given together with the direction of the free 
stream, both the right hand side of the equation and the kernel of the 
integral are known. Hence the equation can be solved for 9(s). In 
principle, the equation may be solved by analytical methods such as 
Neumann series successive approximation. In practical application, 
however, numerical methods are appropriate. Based on this concept, 
the contour of the given airfoil is divided into N segments. The 
integral term in the Integral equation may be written as 

nVi a 

iiiJ <J(«’)^^«r(s,s')ds' 
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where s^ and are the end points of segment j. At this stage, Hess 

and Smith make the approximation that a(s) takes a constant value within 
each segment, This approximation will tend to hi exact when N goes to 
Infinity. With this approximation, the integral becomes 


^ -d 

T 0, 1 ^ An r (s,s* )ds' 

3-1 ' ij 


where Oj is the constant value of a(s) in segment j, In order to carry 
out the integrations analytically and hence simplify the problem further, 
two more approximations are made. One is to approximate the curved 
segment j by a straight line joining the end points of segment j With 
this approximation, each integral may be evaluated for a specified point 
s regardless of the precise shape of the segment j. Since the coordinates 
of the end points of each segment must be known, a convenient choice of 
this point s would be the mid-point of the chord line of each segment 
The second approximation is chat the integral equation (23) is not to 
be applied on the airfoil surface. Instead, it will be applied at the 
mid-point of the chord line of each segment. Therefore, N equations may 
be obtained by applying equation (23) at N of these mid-chord points ■ 

This system of equations may be written in the form of 


N 

E K..a, - F. , i - 1,2 N 

j-1 J ^ 


<24) 


where “ 


U * n. and n, is the unit outward normal vector of chord 

SB 1 1 


1. Is the abbreviation of 


r- 


dn. 


in r(8,8*)ds' 


where s denotes the mid-point of chord i and integration has been carried 


out along the chord j. When i equals J, the value of K^j is 2rr which is 
the coefficient of the first term on the left hand side of equation (23) 
Literally, represents the outward normal velocity at mid-point of 
chord i Induced by a source sheet with uniform strength unity located at 
chord j . When the geometry of the airfoil is given in terms of the co- 
ordinates of N discrete points, the coefficient matrix K^. . can be computed 


The right hand side is then determined once the free stream direction is 
given. Therefore, a solution can be obtained by solving this system 
of N simultaneous linear algebraic equations. This solution is an 
approximate solution to equation (23) evaluated at discrete points. Be- 
cause the integration is carried out along the chord line of each segment 
in obtaining and equation (23) is actually applied at mid-chord points 
which are slightly off the airfoil surface, the segmentation of the airfoil 
surface should be made in such a way that segment size is smaller in the 
high curvature region and larger in the low curvature region of the air- 
foil surface. In other words, more segments are needed near the high 
curvature region in order to obtain better results. This dees not liiq>ly 
that the segment can be made very large when a large portion of the air- 
foil surface is a straight line. The reason is that the variable a(s) 
will not be a constant even for a straight line portion of the airfoil 
surface. Hence the segment size must be kept small even for a straight 
line portion of the airfoil surface order that the first approximation 
made by Hess and Smith can be considered to be a good one. 

With the solution to equation (24) , the velocity tangent to the 
airfoil surface can be computed by combining the tangential conq>onent 




sf-5- 1-^-* .• afev^.: 
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of the free stream velocity and the tangential velocity induced by the 
source sheet. This surface velocity usually does not vanish at the 
trailing edge of the airfoil when the airfoil is generating a lifting 
force. Therefore, an additional set of 0^^ must be superposed to the 
solution of equation (24) in order to insure that the fluid will flow 
off the trailing edge smoothly. This set of is obtained as the 
solution to a system of simultaneous linear equations with the same 
coefficient matri;; appearing in equation (24) but with a different 
right hand side. The right hand side for this purpose is the normal 
velocity at each midchord point induced 1 y a vortex sheet with unit strength 
which has exactly the same location as the approximating source sheet. 

The solution to this system is a source sheet which will induce a 
tangential velocity on the airfoil surface which corresponds to the 
surface velocity distribution due to a unit circulatory flow around 
the airfoil. By combining these two sets of and varying strength 

of the vortex sheet, one is able to satisfy the Kutta condition with 
a circulatory flow of strength T. This F essentially represents the 
circulation generated by the airfoil. With this F, the real tangential 
velocity along the airfoil surface can be computed by combining the 
tangential component of uniform free stream, the tangential velocities 
induced by the two sets of 9^'s. Pressure distribution can then be 
obtained by computing the pressure coefficient from 


C - 1 
P 


'"t 2 
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2. Mar tensen- Jacob ' s Method 

9 

Instead of replacing the airfoil surface by a source «heet, Martensen 
uses a vortex sheet, By requiring that the strength of the vortex sheet 
be identical to the velocity distribution on the surface of an airfoil, 
Martensen was able to show that the interior of the closed vortex sheet 
must have zero velocity everywhere In particu. j.r . the tangential 
velocity at every point on the inner side of the vortex sheet caused by 
the free stream and the vortex sheet should be zero. This result can be 
represented also by a Fredholm integral equation of the second kind 

YfsVan r(s,s’)ds' 

0 

= U cos a + ^ sin a) (25) 

03 ds ds 

where y(s) is the strength of the vortex sheet and a denotes the free 
stream direction This equation almost has exactly the same form as 
equation (23). The corresponding terms possess similar meaning except 
that the tangential velocities are considered in equation (25) 

When solving equation (25), Martensen also replaces the integral 
by a summation, and hence a system of simultaneous equations is to be 
solved. However, when approximating the integral, the first law of the 
mean is used in contrast to the second law of the mean employed by Hess 
and Smith. When 2N control points are distributed along the airfoil 
surface, the resulting system of simultaneous linear equations can be 


written as 




E K^jYj - 2 Nw^ , i - 1,2,...,2N 
where Yj t:he vortex strength at point J and 


“ Ujj[xj^(6)cos a + y^(6)sln a] 

(x^-x^)y^(e) - (y^-y^)x^(S) 


(Xj^-x^)^ + (y^-Xj)^ 


2N 

= - S K 
i=l 
iVj 


i ^ i 


i = j. 


The independent variable has been changed from s to d where d is the 

angular coordinate of the image of a control point when the airfoil is 

transformed into a circle. The control points on the airfoil surface 

are to be distributed in such a way that their images are equally spaced 

on the circle, and hence the interval over which the law of the mean is 

applied is Dots over x and y indicate derivatives with respect 

2N 2N 

to 6 . The identity 2 K = 0 comes from the fact that 


J |— ^n r(s',s)ds' = - tt. 


Equation (27) makes the solution Y(s) of equation (25) non-unique and 
this character is carried over to the system of equations which replaces 
equation (25). Indeed, equation (26) is not a linearly independent 
system and the coefficient matrix has a rank cf only 2N-1. The only 
degree of freed<xn of the system is annihilated by applying the Kutta condition. 
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The solution obtained Is the velocity distribution on the airfoil surface. 
This is an advantage of this method when compared with Hess>Smlth method, 
As previously stated, the direct solution to equation (24) Is the source 
strength, and the surface velocity has to be obtained fron. this source 
sheet. However, Martensen's method does not give good results for thin 
airfoils with thickness/ chord ratios less than 10%. The reason is that 
when the airfoil is thin, the control points on the upper and lower 
surfaces are very close to each other, and the vortices located there 
induce strong tangential velocities on each other. While this Induced 
velocity decays very rapidly for points in the neighborhood of the 
control point, the first law of the mean assumes it to be a constant. 
Therefore, the computed result Yj Is not a good approximation to the true 
solution to equation (25). Jacob^*^ modified Martensen’s method by taking 
the limit of when i approaches j to be the value of . This gives 


K 


ii 


1 x.(0)Vi(0) - y^(e)x.(e.) 

2 ^ 2 ^ 2 

x.(0)^ + y.(9) 


Then the value of for control points directly facing each other,, e.g., 

i " 2N - J + 2 for symmetrical airfoils, is obtained from ^ K. . ■ 0. 

i«l 

This modification improves the results obtained from Martensen's method, 
but it places a restriction on how the control points are to be 
distributed over the airfoil surface. Also, it can cause a lot of in- 
convenience for airfoils such as those obtained by Liebeck and Ormsbee ^ 
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3. Pel lets' Method 

Oellers^^ developed a method to compute the pressure distribution on 
the surface of cascades of airfoils In tandem, the surface of each air- 
foil Is replaced by a vorte>: sheet which must be a streamline. Instead 
of working with induced velocities, a method used by both Hess-Smlth 
and Mar tensen* Jacob, stream functions are employed. The stream function 
for a uniform 'free stream is added to that of the vortex sheet, and the 
sum is set to be a constant on the airfoil surface. This requirement is 
represented by a Fredholm integral equation of the first kind 

1 

^ - u^y(s) - v^x(s) - J Y(s') 

0 

X Xn{sinh2 !l xX s J-x(s»)l ^.^2 l£ y . (s)-y<s')-l ^l/2 

^ere u ■ U cos a» v = U sin a, t is the spacing of the cascade and 

CD 00 00 CD 

ijr is an unknown constant. When a single-element airfoil Is under con- 
sideration, the kernel of the integral becomes in r(8,s') and the equation 
reads 

1 

t • r(8,8')ds' . (28) 

0 

To solve this equation for ^ and y(&)> the Integral Is again replaced by 
a siommatlon using second law of the mean to approximate the Integral. 
Dividing the airfoil surface into N segments and applying equation (28) 
at the mid -chord point of each segment results in a system of 
simultaneous linear equations of the form 
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♦ “ H-yi - +jS]. i 

N 

Z -t - - xx^y^ , 1 - 1,2,...,N (29) 

For an airfoil with the geometry given in terms of coordinates of N points, 
the coefficient matrix may be computed and the right hand side of 
equation (29) is knownonoe the free stream direction is given. However, 
there are Hfl unknown variables? and y^'s while only N equations are 
available. The solution is not unique. This degree of freedom is again 
removed by applying the Kutta condition at the trailing edge. Because a 
vortex sheet is used to represent the airfoil contour, the solution 
is the real tangential velocity on the airfoil surface. This is similar 
to Martensen-Jacob's method but Oellers' method is simpler due to less 
computation involved in obtaining the coefficient matrix. 

4. Comparison of the Three Methods of Computing Pressure 
Distribution on the Surface of a Given Airfoil 

Extensive investigations have been made in order to find various 
properties of the methods presented in the previous subsections and their 
applicability to solving the problem of airfoil design. The results 
are summarized as follows. 

When applied to standard airfoils for which analytical expressions 
of pressure distribution are available, the Hess-Smlth method always 
gives the correct value of circulation generated by the airfoil . This 
indicated that the pressure distribution computed by the Hess-Smith 
method is fairly close to the true value. However, the computed surface 
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velocity is found to be very sensitive to the coordinates of the control 
points on the airfoil surface which are used as the input of this method 
That is to say, the numerical values of the input coordinates have to be 
so accurate that they really do represent a smooth curve A tiny error 
in the input coordinates can produce a wavy behavior of large ao^litude 
in the conqjuted surface velocity This is illustrated in Figure 7. The 
upper half shows the trailing edge of a thin airfoil and the lower half 
shows the pressure distribution computed by the Hess-Smith method. An 
airfoil whose surface is not perfectly smooth will produce wavy pressure 
distribution along the airfoil surface but the amplitude computed by the 
Hess-Smith method is simply too large. This is attributed to the fact 
that this method chooses the velocity normal to the airfoil surface to be 
the variable to work with. Because airfoils are generally operating at 
moderate angles of attack and hence the free stream does not have large 
con^onent normal to the airfoil surface, the consequence of using normal 
velocity as the variable is that a small absolute error introduced in 
computing the unit normal vector n^ by the non-exact input coordinates 
is a large relative error for the right hand side of equation (24), 
Therefore, by failing to represent the airfoil surface exactly with the 
input coordinates, the solution of equation (24) contains errors. These 
errors are especially large near the trailing edge of thin airfoils. 

This character of the Hess-Smith method makes it intractable to perform 
modifications of the geometry graphically when an airfoil design 
problem is to be solved Alternatively, if the modification is to be 
accomplished by an iterative procedure programmed for a con^uter, the 
possibility is eliminated when one recalls that the tangential velocities 
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on th« airfoil surface have to be computed as those Induced by all the 
discrete sources distributed on the approximated airfoil contour. 

As noted by Jacob, Martensen's method does not give accurate results 
for thin airfoils. This drawback has been improved by Jacob by taking the 
steps described in Subsection However, this modification 

restricts the manner in which control points are to be distributed on the 
airfoil surface. As will be seen in the next chapter, all the airfoils 
which generate optimum pressure distributions have a sharp trailing edge 
with which a thin aft-part is Inevitably associated. Therefore, a 
modification to the Martensen-Jacob method has been made during the 
course of this study in order that thin arifoils may be treated and no 
generality of how to distribute control points will be lost. The 
modification is to approximate the integral by using the second law of 
the mean Instead of the first law of the mean. In other words, the 
scheme by which K^j is computed in the Hess -Smith method is applied to 
equation (25). It is found that this modification serves the purpose of 
fulfilling the requirements described above, but the circulation generated 
by the airfoil is smaller than the one coaq>uted by the Hess 'Smith method. 
Further investigation reveals the fact that the influence of local 
curvature on the induced velocity has been ignored when the coefficient 
sutrix is obtained by carrying out the integration alrag the chord line 
of each segmeht instead of along the airfoil surface which is usually 
curved. This approximation did not bother Hess -Smith's results but the 
effects on the solution of equation (25) are far more than would be 
ai^ectad. The reason is a very interesting one: The diagonal almunts 
of the coefficient suitrix of the integral equation (23) have a positive 


sign while those of the integral equation (25) have a negative sign. 

With everything else the same, this difference makes the coefficient 
matrix of equation (24) a diagonally dominated one while that of equation 
(26) is an almost singular one. When one uses the terminology of numerical 
analysis, the former is a we 11 -conditioned matrix and the latter Is an 
ill-conditioned one. Hence by neglecting the local curvature, a more 
severe consequence appears in Martensen- Jacob ' s results Chan In Hess- 
Smith's results. After the curvature effect is taken into consideration, 
the circulation computed by Martensen- Jacob ' s method increases, but it 
is still smaller than the value obtained by Hess-Smlth's mthod. On the 
other hand, because it is the vortex sheet and tangential velocities i^ich 
are considered by Martensen and Jacob, the cor^puted results are not very 
sensitive to the inaccuracies of input coordinates. Small absolute 
errors introduced by non-exact coordinates in computing tangential 
vectors produce small relative errors fo;' the right hand side of equation 
(26) Thus, wavy pressure distributions will be obtained when non-exact 
coordinates are used as input, but this is what should be expected. The 
wavy pressure distribution computed by the Hess-Smith method for a same 
set of coordinates is too exaggerated to be realized in real fluid flow 
Hence the modified Martensen-Jacob method can be used in designing air- 
foils by an iterative procedure dien modifications of geometry are to be 
made graphically. Also, the explicit appearance of surface velocity 
Y(s) in the integral equation maWes it possible that, when applied to 
airfoil design, a syst^naticsl way of modifying the geometry may be 
established and performed on a computer. 
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There are no references which give the pressure distributions along 
the surface of single-element airfoils computed by Oellers* method, but 
the numerical examples which have been worked out during the course of 
this study show that Oellers' method possesses many favorable character- 
istics. First, the confuted surface- velocities are relatively Insensitive 
to Inaccuracies of the input coordinates. This is attributed to the fact 
that the right. hand side of the system of simultaneous equations (29) 
contains only the coordinates of the airfoil while derivatives must be 
conq>uted both in the Hess-Smlth method and the Mar tensen- Jacob method. 
Because computing derivatives numerically always causes a loss of accuracy, 
the results obtained by using Oellers* method are expected to have higher 
precision than the other two methods. Second, the circulation generated 
by the airfoil is found to be almost the same as the one computed by the 
Hess-Smith method. The word 'almost' is used here because all the 
numerical results are approximations to the real solution to the Integral 
equation en^loyed. Therefore, an identity of the results obtained from 
different approximations is almost Impossible to achieve. The pressure 
distributions on the surface of standard airfoils are found to be very 
close to those obtained by analytical methods. This is true whether the 
airfoil Is thick or thin or whether it has a rounded trailing edge or a 
sharp one. Third, the time consumed in conq>uting the velocity distribution 
Is less than that of any of the other methods. This is attributed to 
the slu^llclty of the kernel of integral equation (28). This character- 
istic Is very Important because the computation of surface velocity isust 
be repeated many times during an Iterative procedure of airfoil design. 
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Thest three eherecterlstlcs, plus the fact that surface velocity y(s) 
appears explicitly in the integral equation, make the Oellera' method 
undoubtedly the appropriate tool to be used in airfoil dcslg;n. 

C. A New Method of Two-Dimensional Airtcil design 

As indicated in the previous sections, the problem of multiple* 
element airfoil design can be solved best by an iterative procedure in 
vhich the geometry of a starting airfoil undergoes modifications until 
the surface pressure distribution computed by a reliable method agrees 
with the desired one. The most reliable method ot computing pressure dis* 
trlbution on the surface of an airfoil with given geometry has been found 
to be the one by Oellers while the modification of starting geometry has 
yet to be studied. Since the modification should be made according to 
how the computed velocity distribution differs from the desired one in 
order that the procedure converges to the desired answer, some means 
must be found by which the modified coordinates of the starting airfoil 
can be related to the desired velocity distribution. When the equation 
considered by Oellers in computing the surface velocity of an airfoil Is 
recalled, 

1 

♦ ■ ’ 27 x 1 r(8,s')ds’ (28) 

0 

one can see clearly that this is the appropriate relation between surface 
velocity *>^<1 coordinates (x,y). In order to compute the pressure 
distribution on the surface of an airfoil, the coordinates of the air* 
foil may be normalised with respect to the chord length for the purpose 
of convenience. With this normalisation made and approximately one 
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hundred control points distributed along the airfoil surface i the co- 
efficient matrix obtained by approximating the Integral in equation 

(28) It found to possess an excellent property. That Is: most eleiaeots 

•3 *2 

of are of order 10*'’, some of them are of order 10** and only a few 

are of order 10*^, Therefore, If the equation Is rewritten as 

®T 

■ v^x(s) ■ ♦ + J Y(s')fn r(s,s’)ds’ 

0 

which can be replaced by 


N 

. i = 1.2 N, 

the change of the sunsnatlon term will not be large If each is changed 
by an amount of order one. This nice characteristic of the coefficient 
matrix forms the foundation upon which an Iterative method of designing 
airfoils Is based. 

The method starts from a computation of velocity distribution along 
the surface of an airfoil. The geometry of the airfoil is an arbitrary 
one and It Is given In the form of a set of coordinates. With this set 
of coordinates and a given free stream direction, the surfa^ velocity 
can be obtained by solving the simultaneous linear equations 


N 

E iYj - t • » i"l,2,..,,M (29) 

for Yj'* t* This surface velocity generally does not agree with 
the desired distribution because the coordinates represent merely an 
arbitrary airfoil. Hence the coordinates need to be changed in order to 
obtain the desired velocity distribution. At this tlaie, equation (29) Is 
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satisfied by the coordinates (x,y), the from this set of coordinates, 
Che free stream direction and the computed if- This identity will 

be destroyed when the confuted is replaced by the desired velocity 
If the desired velocities kept there and an attempt is 

made to change the values of x,y in such a way that equation (29) is again 
satisfied, it would mean that an airfoil \diich produces the desired 
velocity distribution has been obtained, and this airfoil is represented 
by these new coordinates. In doing so, u^ and may be kept unchanged, 
is the desired value and it causes no trouble. Now, the equations 
are to be satisfied by varying x, y, y and while is strictly 
determined by (x,y). This puts too many constraints on the effort of 
making attenq»ts to satisfy equation (29) Thus, an alternative approach 
is taken in which the coordinates of the airfoil which produces the 
desired velocity distribution are to be obtained by making several 
modifications to the original airfoil instead of one modification. The 


first step in achieving this is to retain the value of which 

corresponds to the original (x,y). The reason for doing so is that , 

2 

consists of N elements. Each element has a different value, and it is 
impossible to find correct values for all the elements in order to 
satisfy equation (29) . Now, the summation § ^ij^j have two values 
not much different from each other whether the Yj is the desired velocity 
distribution or the one computed for the original airfoil. Therefore, 
it is to be eiqiected that only small changes of x, y and (f are needed 
to restore the identity of equation (29) after Yj^^'s are replaced by 
yj**^*!. Since the goal is to change values of x and y, if is assumed 
to be unchanged and its original value is retained. At this point, the 
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situation Is that the left hand side of equation (29) has been computed 
using the original values of and desired values of Vj's# u^ snd 
are unchanged and the values of x and y are to be found such thst the 
identity can be restored. Since N equations can uniquely determine only 
N variables, either x or y must be forced to take its ortglnel value. 

A study of Figure 8 will help in making this decision. In this figure, 
two airfoils are shown Immersed in a common free stream. These two alr« 
foils have entirely different angles of attack, thickness distributions 
and camber distributions. The only thing they have in common Is that they 
span the same lengths in x direction. If the geometry of an airfoil Is 
represented by a set of control points which are the intersections of 
the airfoil contour and the family of vertical lines shown In the figure, 
the geometry of one airfoil can be obtained simply by moving the control 
points of the other airfoil along these vertical lines. In other words, 
the geometrical characteristics of an airfoil can be completely changed 
by changing the y coordinates of the control points. Therefore, In 
restoring the identity of equation (29), the x values are assumed to have 
the original values. The system, then, contains N y>values to be deter- 
mined from N equations. The solution to this system is easy and the 
y-values obtained will be designated by y^^\ This y^^^ together with 
the unchanged x-values actually represent the control points of an air- 
foil which would produce the desired velocity distribution If the co- 
efficient matrix remained unchanged during the change of y-values. 
Because Is completely determined by (x,y), the new airfoil represented 
by (x,y^^^) will have a different coefficient matrix, Mid hence 

Che desired velocity distribution will not be realized. However, since 
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which hee been found co be smell, the change in is to be expected to 

be smell also. Hence the velocity distribution on the surface of this 

new airfoil is not much different from the desired distribution . 

J 

to this point only one iteration has been completed. The geonetry of 
this new airfoil can be modified again by changing the y^^^ values by an 
amount 
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( 2 ) 

The velocity y' along the airfoil surface which is now represented by 
(x,y^^^) is expected to be closer to y^^^ than y^^^ and y^®^ because the 
difference y^^^ - y|^^ is smaller than y^^^ - y^ in absolute value. There- 
fore, an iterative scheme has been established to obtain the geometry of 
an airfoil with velocity distribution specified along the airfoil surface. 
As will be seen in Chapter V, certain precautions must be taken when 
specifying the velocity distribution in order that a satisfactory airfoil 
may be obtained. For the examples presented in Chapter IV, the con- 
vergence of this iterative process was found to be fairly r^id. For 
clarity, the procedure can be formally stated as follows: 

(i) Qioose one of the standard airfoils as the starting con- 
figuration. Either Joukowski airfoils or HACA airfoils may be used be- 
cause Btanderd fomwlae ere available to compute the coordinates of the 
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control points dtstributsd along the airfoil contour. Hora control 
points arc needed in high curvature regions in order to obtain accurate 
results. The coordinate systen is to be oriented in such a way that the 
chord of the airfoil is parallel to the x-axls, and a free streaai direction 
should be chosen. For a nonaallsed chord length, one hundred control 
points were found to be satisfactory for the exas^les cos^uted in this 
study. 

(ii) With the control points chosen in (i), cos^ute the surface 
velocity by Oellers method. 

(ill) Compare the confuted velocity distribution with the desired 
one. Replace those Yj'a vhlch are not desired by the desired values and 
evaluate 

N 


(iv) 


Obtain the y«coordlnates of the modified airfoil from 


( 30 ) 


where superscripts (m) and (m*l) dmote the sequence of iteration and 
superscript (d) denotes the word 'desired*. 

<v) Compute the velocity distribution along the surface of this 
modified airfoil by Oellers Mthod and return to (ill). 

(vi) The repetition of steps (ill), (iv) and (v) is to be etopped 
when the velocity distribution obtained in (v) is satlsfeetorily close 
to the desired distribution. 
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Th« prlnelplt of this new method of airfoil design Is net eemplieated, 
but seme difficulties are encountered when an exercise of the design 
preeedurs is carried out. This will be seen In Chapter V. 

For the purpose of demonstration and simplicity, the new method of 
airfoil design has been stated for problems of designing slngle*elemeat 
airfoils. The entenslon to the problem of designing multlple-elesMint 
airfoils Is issaedlate. The appearance of several auxiliary elements In the 
flow field. In addition to the laaln elanent. simply adds more Integral 
terms to equation (28), and the value of stream function if will be 
different for different elements. Therefore, equation (28) may be 
rewritten as a system of equations 




P P 


q q P q 4 


p * 1.2. ....n 

where subscripts p and q denote the airfoil elements p and q respectively, 
and the system Is assumed to consist of n elements. If the control points 
of each element are numbered In such an order that they all start from 
the trailing edge and go around the contour of each element In the same 
direction, a system of simultaneous linear equations can be written as 


j-i 


■ ♦- ■ ''-*1 ■ 


1.2...t.If_ 


(31) 


where subscripts i and J denote the control points 1 and j respectively. 

is the total number of cmtrol points and subscript m denotes the 
m-th eleoMnt. There are n imknown t^'s and unknown Yj's to be solved, 
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but only Kj. equations ere evellable. This n degrees of freedom will be 
removed when the Kutte condition is applied at the trailing edge of each 
elasMmt. Therefore, the procedures of designing oMltiple^elenent airfoils 
are identical to the ones for designing single-element airfoils except 
that equation (31) should be used instead of equation (29). 

It is important to point out that although this new method of air* 
foil design is. powerful in the sense that the starting geoisetry does not 
have to produce a velocity distribution which is very close to the desired 
one in order to achieve a converging iterative process, there is a dis- 
advantage in the way the airfoil geometry is modified. As stated in 
step (iv), the geometry of the airfoil is modified by changing the y- 
coordinates of control points. This modification can entirely change the 
character of a single-element airfoil. But, when multiple-element air- 
foils are under conaideration, this method is not capable of making all 
the modifications permitted by the oxlstlng degrees of freedom. This is 
because each element ought tc be able to move freely relative to the 
others during the modification while step (iv) only allows motions in 
the y-dlrection. Therefore, the relative position in x-direction is 
fixed once a starting Is »fven. This disadvantage, as will 

be seen in Chapter V, results in a possibility that an airfoil which 
produces the desired velocity distribution may not be obtained. Never- 
theless, the necessary change of relative positions in x-direction still 
can be s»de by artificial means because the dis^vantage described above 
exists only when the systematic steps (i)--(v) are to be programsMd for 
a computer and the coordinates of crn^rol points are to be changed by 
the awchine. 
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D. Selection of Parameters in DesiRnina Two-Element Optimum 
Airfoils 

As stated In Section II. B, the parameters which specify the optimum 
pressure distribution are the trailing edge velocity and Reynolds number 
RSq. Because a non*zero trailing edge velocity Is necessary to avoid an 
airfoil of Infinite chord length, a convenient choice would be • 1. 
Namely, the velocity at trailing edge Is the same as free stream velocity. 
This choice is. made for all the airfoils generated in this research and 
hence only different values of Re^ need to be considered. As can be seen 
in Stratford's derivations, Re^ specifies the boundary layer character 
Istlcs at Xq and honce the zero skin-friction pressure distribution. In 
the variational problem, with trailing edge velocity fixed at ^ 

Re^ determines the peak velocity, the extent of peak velocity plateau and 
the chord length which corresponds to a specified free stream velocity. 
Hence, to design Sk single-element optimum airfoil, an ^propriate value of 
ROq must be chosen according to the desired free stream velocity and chord 
length. The exact value of liit coefficient is not known until the 
Iterative process Is terminated because the velocity distribution along 
the lower surface does not have a definite specification for the reason 
stated in Section II. C. For single-element airfoils, a given velocity 
distribution on the upper surface can be achieved by different airfoils 
having appropriate coBd>lnatlons of angle of attack, camber distribution 
and thickness distribution. Any one of these airfoils may be considered 
to be the solution to the optimization problem. However, the final decision 
as to which one Is best suited for practical utilization depends on other 
criteria such as structural requlres^t, performance at various angles 
of attack, etc. 
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When two elements are present in the system, things are more 
flexible. To avoid unnecessary complexity, the requirement « 1 is posed 
at the trailing edge of both elements. A value of ROq needs to be 
specified for each element and they depend only on the desired chord length 
of each element. With this information, the pressure distribution can be 
specified, and the iteration can start from an arbitrary configuration. 

In addition to the geometrical characteristics possessed by single '•element 
airfoils such as angle of attack, camber distribution and thickness 
distribution, two'-element airfoils have one more, namely, the relative 
position of two elenents. With the velocity distribution specified 
definitely only on the upper surface of each element, there are many air- 
foils which can be considered as the solution to the optimization problem. 
These airfoils all produce the desired velocity distribution on the upper 
surface of each element, but they will have different angles of attack, 
different thickness and/or camber distributions of each element and 
different relative positions between the two elements. Once again, 
structural aspects and aerodynamic behaviors at various attitudes play 

important roles in determining which airfoil is optimum in the sense 

max 

and utilizable in constructing wing sections. 
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IV. EXAMPLES 

During the development of the new method of airfoil design, a two- 
element airfoil was considered in which one element was placed near the 
trailing edge of the other. This simulates the conventional wing which 
possesses a slotted flap as the high lift device. Values of 5 x 10^ and 
10^ were chosen for Re^ because they give a chord ratio approximately 4. 
After the development of the method was completed, another two-element 
airfoil was generated with Re^ = 10^ and 2 x 10^. With these two sets of 
Reynolds numbers, different relative positions were assigned to both air- 
foils in order to investigate the effects of relative position of the 
two elements. The results are presented in graphic form in Figures 9 to 
18. Both geometry and pressure distribution are shown in the figures. 
Free stream direction is at an angle of 11.25° measured from positive x- 
axis. Each figure shows the result obtained by conqpleting ten iterations 
of the procedure described in Section III.C. The time consumed by an 
IBM cong>uter model 360/75 were approximately 187 seconds and the figures 
were drawn by the Calcomp Plotter. Since the Plotter only draws straight 
lines, symbol x indicates where was computed, and Lagrangian inter- 
polation was employed in plotting the airfoil contours. 

In addition to the two-element airfoils, two single-element airfoils 
were generated. Values of 5 x 10^ and 10^ were chosen for Re^ in order 
to cong>are the geometries «rith two-element airfoils. Sato's method was 
employed and the resulting airfoils were found to possess highly curved 
leading edges. Also, the thickness ratio was found to be less than 9%. 
Shown in Figure 19 is the airfoil for Re^ « 5 x 10^ together with three 




10. Geoaetry and pressure distribution of tno~eleiiient optiana 
airfoil for Re^, « 5 x 10® and Re^, = 10®, C, « 2.32828. 














Figure 16. Geoaetry and pressure distribution of two-eles»nt option 
airfoil for Re„, - 10^ and Re„, - 2 x lO®. c » 2.49528. 






Figure 17. Geometry and preaaure distribution of tvo- 
element optimum airfoil for Re>. ■ 10^ end 
HeQ2 - 2 X 10^, • 2.45446. 



Figure 18. Geometry and pressure distribution of two-element optimum 
airfoil for Re^^ «• 10^ and Re_, - 2 X 10^, C_ = 2.50194. 
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Figure 19. Single>elenent airfoil designed by Sato's 
method and pressure distributions computed 
by three different methods. 
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pressure distributions. These three pressure distributions were obtained 
by three different methods and are presented in the same figure for 
comparison. It can be seen that Oellers' method yields a better result 
than the one obtained from the Hess-Smith method. Since thick airfoils 
with blunted leading edges are desirable for practical utilization, 
modifications in pressure distributions were made, and the new method of 
airfoil design described in Section III.C was employed to obtain the 
coordinates. Ninety control points were used and the time consumed by 
the IBM 360/75 con^uter was 110 seconds. The resulting geometries and 
pressure distributions are presented in Figures 20 and 21. Also shown in 
the figures are the pressure distributions at off-design attitudes, i.e., 
at angles of attack smaller than the designed values, 17.6° for 
ROq « 5 X 10^ and 18.6° for Re^ = 10^. 

Because of the different quantities used In computing Reynolds nuniber, 
the free stream Reynolds number Re may be obtained from the formula 

CO 


Re 

00 


Re 


c 



where all the quantities are available treta Tables 1, 2 and 3. The value 
of c/s^ is approximately 0.925 for the exanqples generated. For wing- flap 
configurations, the main element is taken as the reference. 

Although s^ has been used as the reference length in defining the 
lift coefficient in Chapter II, the chord length c is commonly used. 

Hence the values shown in Figures 20 and 21 are based on the chord 
length and those shown in Figures 9 to 18 are based on the chord length 


of the main element. 






x/c 

Figure 21. S ingle ^element optimum eirfoll for 
Re^ « 10^ and pressure distributions 
at three different angles of attack , 
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V. DISCUSSION AND CONCLUDING REMARKS 

A. Single-Element Optimum Airfoils 

Although the main Interest of this research was to obtain the 
geometries of multiple-element airfoils optimized for maximum lift co- 
efficient, several single-element airfoils were also generated. Two 
single- element airfoils with different Rejmolds numbers are presented in 
Figures 20 and 21 together with their pressure distributions. The 
geometries of these two airfoils are quite similar. This is not surprising 
since the optimum pressure distributions on the upper surface of each air- 
foil do not differ very much (see Figure 22). The modified pressure 
« 

distributions on the lower surface are also similar. The only difference 
between these airfoils is the attitude at which each airfoil is to be 
operated. As can be seen in Tables 1, 2 and 3, larger Reynolds numbers 
always result in higher peak velocity and longer extent of peak velocity 
plateau. Hence lift coefficient increases as Reynolds number increases. 

As a consequence, single-element optimum airfoils with larger Reynolds 
numbers usually operate at higher angles of attack The geometry is 
obtained first by Sato's method. The results all possess leading edges 
of high curvature and the thickness ratios are smaller than 9%. This 
highly curved leading edge is the consequence of the abrupt pressure drop 
at the front stagnation point as demanded by optimization of pressure 
distribution, but it is not suitable for operating the airfoil at other 
angles of attack. When the airfoil is to be operated at angles of attack 
smaller than the design value, the front stagnation point moves along the 
airfoil surface in a clockwise direction. This leads to a small region 
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of high velocity on the lower surface and consequently a flow separation 
when the velocity is decreasing toward the trailing edge. Although the 
flow separation on the lower surface is probably not as extensive as the one 
on the upper surface, the precise effects are not known. The curvature . 
near the leading edge was reduced by modifying the geometry obtained from 
Sato's method in an arbitrary manner. This modification not only changed 
the abrupt pressure drop to a gradual one but also destroyed the optimum 
pressure distribution which is composed of a plateau and a Stratford zero 
skin friction pressure distribution What needs to be done is to modify 
this modified airfoil in such a way that the gradual pressure drop, as 
well as the plateau and Stratford's dittribuuion, can be obtained. This 
modification is accomplished by employint'* the new method of airfoil 
design developed in Section III.C, As shewn in Figures 20 and 21, each 
airfoil has a leading edge of small curvature, a moderate thickness (127.) 
and a not too thin aft part. The pressure distribution on the upper 
surface does have a plateau and the recovery part does follow the 
Stratford distribution. But, the pressure drop at the front stagnation 
point is gradual and the trailing edge velocity is not - 1. These 
facts make the airfoil a non-optlmized one. However, in view of the 
conflict between the optimization requirement and other criteria such as 
the aerodynamic performance at various angles of attack and structural 
requirement, the compromise made is an appropriate one. 

The behavior of these airfoils at lower angles of ar.tack was also 
examined. The pressure distributions are shown in Figure 20 and 21. 

They all possess an almost constant free stream pressure on the lower 
surface. The acceleration on the upper surface is gradual and the 
position of peak velocity is fairly close to the mid>chord. 
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B. Two-Element Optimum Airfoils 

The examples presented in Figures 9 to 18 are two-element optimum 
airfoils with free stream Reynolds numbers in the order of 10^ to 10^ 

Because the intention is to design high lift devices, such as a retractable 
trailing edge flap, chord ratios of approximately 4 are used. It can be 
seen that geometries of the main elements differ from those of the single- 
element airfoils with the same value of RSq. This difference occurs be- 
cause the optimum pressure distribution is the same for all airfoils of 
any number of elements while for multiple-element airfoils the pressure 
distribution of each element is not entirely determined by its own geometry. 
The geometries of ocher elements an^ the relative position of each element 
are also involved. These examples have indicated that the relative 
position of two elements plays an important role in determining the 
geometry of each element. VJhen two elements are in a conventional wing- 
flap configuration as shown in Figure 9, the influence of the flap on the 
main element is to induce an increase of circulation which is composed of 
a velocity increase on the upper surface and a velocity decrease on the 
lower surface of the main element. However, in order to produce the same 
optimized pressure distribution, the single-element airfoil is inclined 
more to the free stream than the main element of a two-element airfoil is 
On the other hand, the influence of the main element on the flap is to 
decrease the circulation about the flap. This is indicated by a decrease 
of velocity on the upper surface of the flap. Hence the inclination of 
the flap with respect to the free stream is more than that of a single- 
element airfoil which produces the same pressure distribution. When the 
distance between two elements is increased, the interaction becomes weaker . 


The consequence is thst in addition to the slight geonetrical change on 
each element the main element increases its angle of attack, and the flap 
decreases its angle of attack. Because the high peak velocity plateau on 
the upper surface of the flap induces a strong flow field in the neighbor* 
hood of the flap's leading edge, low speed near the trailing edge of the 
main element as required by the optimization can not be realized if the 
distance between the two element is too small (see Figures 9 and 14). 
Therefore, from an aerodynamic point of view, some minimum distance must 
be maintained between the two elements in order that the optimum pressure 
distribution can be realized by geometrically realistic airfoils. Placing 
the flap directly beneath the aft part of the main element should be 
avoided. Under that circumstance, the lower surface of the main element 
is approximately parallel to the flap and the peak velocity plateau on the 
flap has a strong Induced velocity. The aft part of the main element 
\.ill contribute negative lift force unless the distance between two elements 
is made large (see Figures 12 and 17). It may be noticed that the require- 
ment q^ ■ 1 at the trailing edge oi the main element is not quite 
satisfied in most cases. This is attributed to the fact that during the 
modification of the airfoil geometries, the relative position of two 
elements in the x direction is kept unchanged. This restricts the freedom 
of modification because the relative position in the x direction is fixed 
at its starting value. Therefore, several starting configurations are 
used and only the one which gives the result of q^ * \ on main element 
should be retained. An alternative solution would be to change the entire 
velocity distribution according to the real value of q^. But, since 
optimum velocity distribution has been reached over a large portion of 
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the airfoil surface and a new value of would result in a new value of 
s^, the gain of changing the velocity distribution is overshadowed by 
the trouble it causes. With regard to the starting configuration, the 
single-element optimum airfoils obtained by a combined use of Sato's 
method and the new method developed in Section III.C is chosen, although 
any geometry may be used. The reason for doing so is two- fold. First, 
all of the single-element optimum airfoils have moderate thicknesses md 
reasonably blunt leading edges. Second, in this way the time consumed by 
the computer in modifying the geometries will not be as long as the time 
which would be consumed by starting from an arbitrary geometry. 


C. Remarks About the Design Procedure 

When the velocity distribution on the surface of a given airfoil is 
to be conq)uted by Oellers' method, two sets of points are used. One is 
the given control points which represent the airfoil. These points are 
on Che airfoil surface. The other set is the set of mid-chord points 
which are slightly removed from the airfoil surface. The integral 
equation (28) is applied at the mid-chord points and the control points 
serve the purpose of constructing the coefficient matrix K^j. Therefore, 
when the geometry of an airfoil is under modification in order to obtain 
the optiimim pressure distribution, the values of y conq>uted by equation 
(30) are the ones for the mid-chord points. In order to compute the 
velocity distribution on the surface of the modified airfoil, a new set 
of control points must be obtained. Two ways of determining the co- 
ordinates of control points from the mid -chord points were considered in 
the early stage of this research but they did not give satisfactory 
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resulci. Okie of these was to take the middle point of the straight line 

which connects two adjacent mid-chord points to be the control point. 

The other was to form a polygon such that the new mid-chord points were 

the middle points of each side and the comers of this polygon were taken 

to be the new control points. Because the x coordinates of the mid-chord 

points were kept unchanged during the modification, these two methods gave 

bad results after a certain number of Iterations. It was felt that the 

means by which mid-chord points were obtained as Initiated by Hess and 

Smith is not appropriate for the purpose of modifying geometries. Between 

the set of control points and the set of mid-chord points, one Is permitted 

to go only from the former to the latter, but not the reverse. To 

liberalize this restriction, interpolation would be a desirable tool. 

« 

In other words, when the control points which represent an airfoil were 
given, an Interpolation was made to give the coordinates of 'mid-chord 
points' where equation (28) Is to be applied. After the new y coordinates 
of 'mid-chord points' were obtained from equation (30), the new control 
points were determined by Interpolation. A Lagranglan four point Inter- 
polation was used and the results were very good. Because none of the 
interpolation methods gave good results vdien the slope of the curve 
formed by the given points was large, certain rotations of coordinate 
axes were necessary when doing interpolation near the leading edge of an 
airfoil. 

During the modification of airfoil getxnetrles, certain values of Yj 
are to be used In equation (30) as the desired velocities. Usually, the 
starting configuration is chosen so that the velocity distribution on 
the lower surface Is approxlsMtely a linear acceleration from leading 
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edge to trailing edge. Hence, those Yj'® which are computed from step 
(ii) of Section III.C as the velocities on the lower surface may be 
retained. Only those on the upper surface need to be changed. In 
generating the examples shown in Figures 20 and 21, it is found that 
too many restrictions on Yj'® result in either an unacceptable geometry 
or a non -converging iterative process. The most severe situation takes 
place when attempts are made in order to reach the exact optimum velocity 
distribution near the trailing edge. The action taken to overcome this 
difficulty is not to pay much attention to the trailing edge. Because 
the viscous effects have not been taken into consideration in the design 
process, the airfoils obtained will not generate exactly the desired 
pressure distribution. This deviation of real pressure distribution from 
the desired one is especially large near the trailing edge. Hence there 
is really no need to attempt to reach the exact pressure distribution near 
the trailing edge. In the case of two-element airfoils, this problem 
seems to be less severe because geometries which correspond to the desired 
velocity distributions can always be obtained by appropriate choice of 
relative position of the elements. 
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